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List of Symbols 


a axial coordinate of duct trailing edge in moving frame 

b axial coordinate of duct leading edge in moving frame 

c ambient sound speed 

k nondimensional characteristic wave number 

L d ratio of duct length to diameter 

m harmonic number 

M = V/c flight Mach number 

M tip - r D Q./c tip Mach number (based on duct radius) 

N B number of point dipoles 

N l number of liner segments 

N obs number of observers for TBIEM3D output 

P total acoustic pressure 

P m m -th coefficient of total pressure 

P i incident acoustic pressure 

P” m -th coefficient of incident pressure 

P t scattered acoustic pressure 

P™ m -th coefficient of scattered pressure 

( r,y/,Z ) cylindrical coordinates in frame of reference attached to duct 

r D duct radius 

r 0 radial coordinate of spinning point dipoles 

t time 

T thrust from fan 
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V 


duct speed 


a(Z) 

M 

e 

K 

Po 

a(z) 

l(Z) 


axial locations of liner segments 

= p 0 c[£(Z) -i<j(Z)] segmented, specific acoustic admittance on interior duct wall 
piecewise specific acoustic admittances 

= 'J I - M 2 compressibility (stretching) parameter 
= k/(} nondimensional stretched characteristic wave number 
ambient density 

segmented, acoustic susceptance on interior duct wall 
segmented, acoustic conductance on interior duct wall 


(O frequency of oscillations (radians/second) for nonspinning dipoles 

£2 shaft speed (radians/second) for spinning dipoles 
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Introduction 


This document describes the ducted fan noise prediction computer program TBIEM3D (Thin duct, 
Boundary Integral Equation Method, 3 Dimensional). The scattering of fan generated noise by a finite 
length, infinitesimally thin circular cylinder in a uniform flow field is considered. The program, based on a 
boundary integral equation method (BIEM), calculates circumferential modal coefficients of the acoustic 
pressure at user specified field locations. TBIEM3D features include versatility, rapid calculations, and 
ease of use. Theoretical and computational details can be found in references 1 -4. 

In a frame of reference attached to the duct, the fan generates spinning acoustic modes. The thrust 
component of fan loading noise is approximated by a collection of spinning point thrust dipoles. A 
precise mathematical representation for the acoustic field due to this configuration has been implemented. 

TBIEM3D employs cylindrical coordinates in a frame of reference attached to the engine (figures 1 - 
3). The coordinate origin is at the center of the fan disc. The fan and duct translate in the +Z (axial) 
direction with uniform speed V . N B equally spaced blades comprise the fan. The shaft rotates with 

speed Q (figure 2). 

The total acoustic pressure in the sound field is split into known incident and unknown scattered 
parts: 

P'(r,\f/,Z,t) — P'(r,\(/,Z,t)+ P'[r,\f/,Z,t) . (1) 

Assuming linear conditions, all dependent acoustic variables can be expressed as superpositions of 

spinning modes. For example, the scattered pressure has the form 

. (2a) 

m=— < oo 

Incident and total acoustic pressures are written similarly. Modal amplitudes are calculated term by term. 
The TBIEM3D code must be run separately for each desired mode. 

For some applications, the propagation of a plane wave through the duct is of interest. The spinning 
point source model is incapable of producing the plane wave. If nonspinning sources are used, then all 
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components of the acoustic field have the form 


P , (r,\f/,Z,t) = e u * f j P m (r,Z)e- imN “ v . 

m=~oo 

The modal function P°{r,Z ) contains plane wave information (see Example 4). 


(2b) 


In TBIEM3D, the point thrust dipole sources may be spinning or nonspinning. The nonspinning 
source configuration could be used to model stators and is useful in examining some aspects of duct 
radiation and propagation. 

The duct exterior is hard and the interior may be hard or lined. The duct liner is modeled by an 
axisymmetric, locally reactive, segmented liner with user specified admittances. The definition of specific 
acoustic admittance used by TBIEM3D 

«(Z) = Po\^7) ~ i<y(Z)] (3) 

is consistent with the time factor e +mNl,Lh in (2). Regions of the duct interior near the leading and trailing 
edges are assumed hard (figure 3). Any interior wall segment may also be rigid. These comments are 
summarized by the equation 


«(Z)= 


(4) 


0 Ze[a,Z,]u[z„ L „,b] 

cx i Z e (Zy,Z y+/ ) j = 1,...N L 
BIEM methodology is a three step process: Step 1) A Helmholtz boundary value problem (BVP) for 
the modal coefficients in (2) is derived. Step 2) Using layered Helmholtz potentials, the BVP is 
converted to a boundary integral equation formulation that features a set of hypersingular integral 
equations for the unknown Helmholtz layers. Step 3) The integral equations are solved and the acoustic 
field calculated from the Helmholtz potential representation. 

The TBIEM3D code is written in the FORTRAN programming language and employs IMSL 
mathematical library routines. TBIEM3D should be implementable on any computer that can 
accommodate FORTRAN and IMSL. Some code modification may be required. For minimally adequate 
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computational performance, a Pentium 133 processor (or equivalent) with 32 megabytes of RAM is 
recommended. 

TBIEM3D input is relatively simple. Geometric, kinematic, and liner parameters are required. If a 
source description other than the one described above is desired, then the user must supply FORTRAN 
subroutines for the calculation of the incident field and its radial derivative. Output from TBIEM3D 
consists of the modal coefficients of the complex pressure components [see equations (1-2)] at user 
specified field points. Postprocessing of results is left to the user. 

The key feature of TBIEM3D is the ability to compute any portion of the sound field without the 
need to calculate the entire field. Competing methods such as finite differences and finite elements lack 
this property. Other positive attributes include reduced consumption of computational resources, 
enhanced numerical accuracy, versatility, coupling of radiation and propagation both forward and aft, and 
validity over a wide range of frequencies. Consequently, the TBIEM3D code is well suited for 
parametric calculations. Many engineering studies of interest can be handled by TBIEM3D. 

Questions, comments, and requests for discussions should be addressed to mhd314@aol.com 
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Limitations 


1) At present, the TBIEM3D code can treat “small” Mach number inflow. Results obtained for 
M > 0.4 may be questionable. TBIEM3D with no inflow restrictions will be made available when 
complete. 

2) For large values of K , TBIEM3D computational time and storage requirements can increase 
considerably. Therefore, at typical fan operating conditions, it is recommended that the user calculate 
a maximum of three circumferential modes. Efforts are underway to improve TBIEM3D 
performance for high frequencies. 

3) It is well known from the theory of wave propagation in an infinite, hard walled duct that resonance 
occurs at certain discrete frequencies. At these eigenfrequencies, the infinite duct problem is 
unsolvable. Theoretically, the finite, hard walled interior duct is solvable at all frequencies. Ill- 
conditioning in the TBIEM3D numerical system, however, is experienced at and near the infinite duct 
eigenfrequencies. TBIEM3D results at these eigenfrequencies show evidence of resonance but 
appear plausible. The numerical correctness of TBIEM3D at resonance has not been established. 
Therefore users should examine TBIEM3D results carefully when the hard wall interior option is 
activated. 

4) For some applications, it may be convenient to place the sources outside the duct. This is easily 
achieved with TBIEM3D. The user must have either a>0 , b <0 , and/or r 0 > r D . 

5) Since the duct is approximated by an infinitely thin cylinder, the acoustic pressure is discontinuous 
across the duct surface. Consequently, evaluation of the acoustic pressure on the duct wall is 
ambiguous. It is recommended that if the pressure on the interior duct surface is required, then the 
user should place the observer a small distance off the duct toward the interior. 
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6) If the nonspinning point dipole option is activated ( kspin = 0), then several TBIEM3D input 

parameters assume different meanings: 1 ) N B = number of nonspinning dipoles; 2) T = arbitrary 

source strength in Newtons which need not correspond to the fan thrust; 3) RPM = {30(0)/ K (shaft 
RPM is meaningless in this context) 

7) The spinning point thrust dipole configuration yields a simplified approximation for the thrust 
component of the fan loading noise. As a result, the user may need to experiment with the radial 
location of the source and the source strength to obtain meaningful quantitative results. 

8) Other sources of fan noise such as thickness and the drag component of loading are not included in 
this version of TBIEM3D. Future versions will include point and/or line source modeling of these 
phenomena. The capability to predict rotor/stator interaction is also being considered for future 
release. 
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Operating Instructions 

TBIEM3D operating parameters consist of a one line identifier, output file name and path, and 
physical parameters. The code generates one output file containing values of program parameters and the 
complex modal coefficients of incident, scattered, and total pressure at user specified field points. The 
output file is associated with logical unit 9. Access of unit 9 elsewhere in the calling program can lead to 
errors and should be avoided. COMMON statements in TBIEM3D should be examined to avoid 
conflicts with the user program. 

To activate TBIEM3D, the user’s calling program must have the FORTRAN statement 


CALL TBIEM3D( kspiti, idem , outfile , m, N B , RPM , r D , a , b , r 0 , T , c , p 0 , V , (5) 

■ M2: ■ {«<£ ■ ’fr}",:- > 


{ 1 N 

x jj J=l in (5) denotes a one dimensional array of length N . Variables in the argument list 


are described in table 1 . SI units are required for dimensional variables. 

The TBIEM3D output file contains the case identifier and program parameters followed by N OBS 
formatted lines containing the dimensional (pascals) complex modal coefficients of incident, scattered, 
and total acoustic pressure. For each observer point, TBIEM3D writes the observer coordinates and 
pressure components according to the following FORTRAN statements: 


WRITE(9,600) Z , r , /?e[/f] , , Re[P s m ] , In{P s m ] , Re[P m ] , In{p m ] 

600 FORMAT(8E11.4) . 
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Examples 

Three examples are presented in this section to demonstrate TBIEM3D features and usage with 
spinning point dipole sources. Kinematic parameters were chosen to simulate actual ducted fan engine 
operating conditions ( M F = 0.40, M TIP — 1.22, N B = 20, h D = 0.50 ). The examples differ in the 
acoustic treatment on the duct interior. Admittances for the three cases are given below. All calculations 
were performed on a Pentium 133 laptop computer with 32 megabytes of RAM. Graphical results 
displayed here are not part of TBIEM3D. 

For each of the three cases, continuous two dimensional portions of the sound field are computed. 
Acoustic pressure and sound pressure level contours for the first harmonic (m = 1 ) are plotted in figures 
5-7. The specific acoustic admittances used for the calculations do not necessarily correspond to actual 
conditions but were chosen for demonstrative purposes. Figure 4 contains the FORTRAN calling 
program that generated the results for figures 5-7. 

Example 1 Hard inlet and hard exhaust. See figure 5. Four minutes computational time required for 
15000 field points. 

Admittance: tt(Z) = 0 Z e \a,b\ 

Example 2 Lined inlet (one segment) and hard exhaust. See figure 6. Eight minutes computational time 
required for 15000 field points. 


Admittance: «(Z) 


1 — i Ze (0,0.4) 
0 elsewhere 


Example 3 Lined inlet (one segment) and lined exhaust (one segment). See figure 7. Eight minutes 
computational time required for 1 5000 field points. 


Admittance: 


\l-i 


a(Z) = 


0.5 

0 


Z e (0,0.4) 
Z e (-0.4,0) 
elsewhere 
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In figure 8, comparisons between the three cases in both the nearfield and farfield are displayed. 
Sound pressure levels for the first modal coefficient on an arc of 200 field points about the duct center 
are calculated. The radius of the arc is ten meters for the farfield example and one meter for the nearfield 
example. The arc extends from the forward duct axis to the aft duct axis. The results are plotted to 
show the effects of liner treatment. Calculations required approximately one minute. 

Example 4 was designed to simulate a no inflow, open ended impedance tube experiment. Twenty 
nonspinning point dipoles were placed near the left end of a three meter long tube of radius 254 cm. A 
one meter long liner specimen with specific acoustic admittance 22 + 02i was placed in the center of the 
tube. The remainder of the tube was unlined. The source frequency was 900 Hz and the zero-th 
harmonic was examined. For these parameters only the plane wave propagates in the unlined portions of 
the tube. 

The complex acoustic pressure for a line of observers at r = 2.0 cm extending along the length of the 
tube was calculated. In figure 9, the amplitude (in decibels) and phase of the complex pressure is plotted 
as a function of distance along the tube. There is approximately a 40 dB drop in sound pressure level 
along the lined section. The complex wave structure in the right end of the tube is due to reflection from 
the open end. 
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Variable 

Description 

Comments 

kspin 

I Integer; source parameter; kspin = 1 implies 
spinning sources; kspin = 0 implies 

1 nonspinning sources 

rr 

ident 

I Character*80 variable; case identifier 

80 ASCII characters maximum 

outfile 

Character*80 variable; output file name 

80 ASCII characters maximum; file path may 
be included 

m 

Integer; Harmonic number 

See (2) and limitation 2 


Integer; Number of fan blades 

See figures 1-3; N B > 0 ; see limitation 6 

RPM 

Real; shaft speed (revolutions per minute) 

£2 = (tt RPM )/30 ; see limitation 6 

r D 

Real; duct radius (meters) 

See figures 1 -3 

a 

Real; axial coordinate of duct trailing edge 
(meters) 

See figures 1 -2 

b 

Real; axial coordinate of duct leading edge 
(meters) 

See figures 1-2 

r o 

Real; radial location of spinning dipoles 
(meters) 

See figures 2-3; see limitation 7 

T 

Real; thrust from fan (kilonewtons) 

see limitation 6 

c 

Real; Ambient sound speed (meters per 
second) 


Po 

Real; Ambient density (kilograms per cubic 
meter) 


V 

Real; Engine speed (meters per second) 

V jc < 1 ; see limitation 1 

N L 

Integer; number of liner segments 

For hard wall interior set N L = 0 in which 
case limitation 3) may apply 

Kir 

Real; axial coordinates of liner segments 
(meters) 

See (3-4) and figure 3; if N L = 0 , then omit; 
a <Z t < Z i+1 < b i = 


Complex; acoustic admittances (non 

dimensional) for segmented liner 

See (3-4) and figure 3; if N L = 0, then omit; 
some segments may be hard, i.e., a i - 0 

N 

iy OBS 

Integer; number of observers for output 

0 < Nobs - '■> Large values of N OBS can 

lead to excessive computational time 

Kir 

Real; axial coordinates of observer points 
(meters) | 

See limitation 5 

bir 

Real; radial coordinates of observer points II 
(meters) | 



Table 1: TBIEM3D Input Parameters 
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Figure 1 : Duct Geometry and Coordinate Definitions 
Cylindrical Frame Fixed to Duct 



0 Plane 
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Figure 2: Source Plane Geometry 
View from Inlet into Duct 


Duct Exterior Wall 
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Figure 3: Duct and Liner Geometry (Side View) 
Two Piece Segmented Liner 


c 

C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 
c Calling Program for Examples 1-3 (Figures 5-7) in User 
c Document - This is a sample calling program and is not part of 
c TBIEM3D 

C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 

c 

program my B I EM 

dimension ZOBS(IOOOOO) , rOBS ( 100000 ) ,Zliner(100) 
complex admit (100) 
character*80 ident , outf ile, outfx 
c 



c Set up TBIEM3D geometric and kinematic input parameters 



c 

ident = 'TBIEM3D Figures 5-7 (see User Document) ' 
outfile = ’c:\xxxl.txt’ 
outfx = 'c:\figs5-7.dat' 
pref = 2 . *10 . ** ( -5 . ) 
nharm = 1 
nblades = 20 
rpm = 3500. 
radiusd = 1.0 
zte - -0.50 
zle = 0.50 
radiusO = 0.90 
thrustc = 1. 
sposnd = 300. 
density = 0.4 
V = 120. 
kspin - 1 
c 



c Input observer points for 2-D field calculations 

c 

c 

nz = 150 
nr = 100 
NOBS = nz*nr 
zmin = -2 . 
zmax = 2 . 

dz = ( zmax-zmin) /float (nz-1) 
rmin = 0 . 
rmax = 3 . 
i f (nr . ne . 1 ) then 

dr = (rmax-rmin) /float (nr-1 ) 
else 
dr = 0. 
endif 
kount = 0 
do i = 1 , nz 

z = zmin+ ( i-1 . ) *dz 
do j = l,nr 

kount = kount +1 
r = rmin+ ( j-l . ) *dr 
ZOBS ( kount) = z 
rOBS (kount) = r 
enddo 
enddo 


Figure 4: Calling Program for TBIEM3D Examples 1-3 
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c 

c 

c Example 1: Hard Wall Interior 

c 

c 

NL = 0 

call TBIEM3D (kspin, ident , outf ile , nharm, nblades , rpm, radiusd, 

1 zte, zle , radius 0 , thrust c , sposnd, density , V, NL, 

2 Zliner , admit , NOBS , ZOBS , rOBS) 

c 

£******************★*********★************************************★** 
c Begin output for user document graphics 
£****★************★************************★************************* 

c 

open (unit=10 , f ile=outfx, status^ ' unknown ' ) 
wri te ( 10 , * ) ' zone t = "1", i = ' f nr, ' , j = ' ,nz 
rewind ( 9 ) 
do j = 1,19 
read (9, *) 
enddo 

do iobs = 1 , NOBS 

read (9,600)zz,rr, res , ais , rei , aii , ret , ait 
pmag = sqrt ( ret * *2+ait * *2 + 1 . e-12 ) 
spl = 20 . *alogl0 (pmag/pref ) 
write ( 10 , * ) zz , rr , spl , ret 
enddo 
c 

£***********************************★******************************** 
c End postprocessing for Example 1 

c******************************************************************** 

c 

c 

c 

c Example 2 : Lined inlet with one segment 

c 

c 

NL = 1 

Zliner ( 1 ) = 0.0 
Zliner ( 2 ) = 0.4 
admit (1) = cmplx ( 1 . , -1 . ) 

rewind ( 9 ) 

call TBIEM3D (kspin, ident , outf ile , nharm, nblades , rpm, radiusd, 

1 zte, zle, radiusO, thrustc, sposnd, density , V, NL, 

2 Zliner , admit , NOBS , ZOBS , rOBS ) 
c 

£********************************** + **+**************************+*** 
c Begin output for user document graphics 
c ******************************************************************** 

c 

write ( 10 ,*) 1 zone t = " 2" , i = ' ,nr,', j = ' , nz 
rewind ( 9 ) 
do j = 1,19 
read (9, *) 
enddo 

do iobs = 1 , NOBS 

read (9, 600) zz,rr, res , ais , rei , aii , ret , ait 

pmag = sqrt (ret**2+ait**2+l . e-12 ) 

spl = 20 . *alogl0 (pmag/pref ) 

write ( 10 , * ) zz , rr , spl , ret 

enddo 

Figure 4 (Continued): Calling Program for TBIEM3D Examples 1-3 
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c 

C*****************************************************************^^ 

c End postprocessing for Example 2 

c* ****************************************************************** * 

c 

c 



c Example 3: Lined exhaust with one segment and lined inlet 
c with one segment 



c 

NL = 2 

Zliner(l) = -0.4 
Zliner(2) = 0.0 
Zliner(3) = 0.4 
admit (1) = cmplx(.5,0.) 
admit (2) = cmplx ( 1 . , - 1 . ) 
rewind ( 9 ) 

call TBIEM3D ( kspin, ident , out file, nharm, nblades , rpm, radiusd, 

1 zte, zle, radiusO, thrustc, sposnd, density, V, NL, 

2 Zliner , admit , NOBS, ZOBS, rOBS) 
c 

*************************** ************************************ 
c Begin output for user document graphics 
c ******************************************************************** 

c 

write (10, *)' zone t = "3”, i = ’,nr,', j = 1 , nz 
rewind ( 9 ) 
do j = 1,19 
read ( 9 , * ) 
enddo 

do iobs - 1 , NOBS 

read (9,600)zz,rr,res,ais,rei,aii, ret , ait 
pmag = sqrt (ret* *2+ait**2+l . e-12 ) 
spl = 20 . *alogl 0 (pmag/pref ) 
write { 10 , * ) zz , rr , spl , ret 
enddo 
c 

c******************************************************************** 

c End postprocessing for Example 3 

c ******************************************************************** 


c 

stop 

600 

format ( 8ell . 4 ) 
end 


Figure 4 (Continued): Calling Program for TBIEM3D Examples 1-3 


20 




C/5 

3 

3 

XJ 

X 

W 

*3 

sc 

TD 

3 


JD 

C 

HH 

*3 

5a 

sc 


4t= 

JO 

£ 

3 

X 

w 

Q 

m 

£ 

W 



«n H 

C/5 

o 

<D 

O 

M S 

l-« 

3 

j Q °^ 

O 

^ o 

00 

ii 

JL> 

3 

<N = 

CN 1-1 

Oh 

_h‘ 

s 

" £ 

w 

g 

o 

Cu 

J°. 

ii 

W) 

Q 

.£ 

O n 

*£ 

Tf 

3 

® o 

Oh 

II <N 

00 

uu II 

l*H 11 


CO 


tn 

<D 

Uh 

3 

tt) 

tC 


21 









Re[P] (Pascals) , SPL (dB, re 20\i Pa) 


CO 






22 


JsPI O' I = X U16'0= J UI0 [= U 4 c 

£'o= u i zi'i = '" 1 n ovo 






a> 

3 

3 

43 

X 

W 

<U 

c 


o 


^ S 

o 
g 2 

HH O 

<D <U 

C 7? 
•-H O 

hJ ^ 


m w 
*= e 

^ O 
Oh Oh 

£ oo 


in 

o 

n £ 




<N 

<N 


On 


II g 

0-t _v 

p9 


3 

X 

W 


Q 

t-H 


e 

c 


Q *& 

ro C/3 


O 

05 o 

II CN 


w 


i> 

<u 

Uh 

3 

Ctf) 

• *«H 

tin 


II 

£ 


23 







co c 
X □ 


C/) 

3 

co 

sz 

X 

LU 


T5 

CD 

C 


~o 

c 

CO 

a> 

c 

"O 

0 

c 


24 


Figure 8: Effect of Lined Interior on Radiated Sound 
Spinning Point Dipole Sources 
M h = 0.4 M tip =1.22 L d = 0.5 

n = l N n = 20 r n =1.0mr =0.9m T=1.0kN 
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Figure 9: TBIEM3D Example #4 - Impedance Tube Simulation 
= 0 r= 0.0254 m M F = 0.0 a = 2.2 + 0.2 i co = 900Hz 



TBIEM3D History 


1) May 15, 1997 : Version 1 released 

2) July 18. 1997 : Version 1 . 1 released; Version 1 errors corrected; nonspinning point dipole sources 
added to accommodate zero-th circumferential mode; User manual modified; Example 4 (Figure 9) 
added; Version 1 no longer supported 
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